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In this paper we study a class of contour dynamics equations depending 
on a parameter a for which has been provided numerical evidence of self- 
similar collapse |2]. This family of equations connect the vortex patch 
problem of the 2D Euler equations (limit a — > 0) with the surface quasi- 
geostrophic equation {a — 1). We explore numerically the evolution of 
these equations, but in order to transform the finite time blowup in an 
asymptotic behaviour, the equations are expressed in self-similar variables. 
We discuss the role of exact self-similar solutions that act as geometrical 
separatrices of what numerically is distinguished as collapsing and non 
collapsing initial conditions. 

Mathematics Subject Classification: 65M99, 35Q35 

1 Introduction 

A classical open problem in mathematical fluid mechanics is whether 3D Eu- 
ler equations may develop singularities in finite time 1 . One of the scenarios 
proposed in the past for the formation of singularities in Euler equations is 
the vortex patch problem. A vortex patch consist of a 2D simply connected 
and bounded region of constant vorticity which is a weak solution of the 2D 
Euler equation. The work of Chemin [5j and Bcrtozzi and Constantin [3] rig- 
orously proved global existence of regular solutions for this case, and therefore 
in the context of this equation, singularities cannot appear. A family of con- 
tour dynamics equations connect this scenario to the evolution of patches in the 
surface quasi-geostrophic equation. This problem, the a-patches problem, has 
been analysed in [2] by Cordoba et al. It is described by the equation. 



Here the sum extends up to P, the number of patches. This is a family of 
equations depending on the parameter a. In the limit a — > the equation 
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approach the vortex patch problem of the 2D Euler equation and the limit a 1 
approaches the surface quasi-geostrophic equation. The later limit has been 
studied mathematically for more than a decade mainly because of its similarity 
to the 3D Euler equations [B] and because being a 2D model it is more tractable 
than the full 3D Euler system. 

The solutions of the a-patches problem have been numerically explored in the 
range 0.5 < a < 1 (see the work by Cordoba et al. [2]) and it has been reported 
evidence of self-similar collapse. In this article we provide a detailed descrip- 
tion of the numerical method used there to perform those simulations and also 
explain how it is applied to the study of the system expressed in self-similar vari- 
ables. The advantage of using rescaled variables is that the finite time collapse 
observed by Cordoba et al, becomes an asymptotic limit in the pseudo-time 
r, and thus the kind of evolution may be more accurately described. We have 
found a family of exact self-similar solutions to the rescaled system. Numerically 
is found that these solutions act as geometrical separatrices of collapsing and 
non collapsing initial conditions. This article describes these findings as follows. 
Section 2 describes the problem. In particular the equations under study are 
introduced and the derivation of the self-similar equations is explained. Section 
3 reports on the numerical methods employed in the simulations. Section 4 
accounts on the results. Finally section 5 presents the conclusions. 



2 The equations 

For the 2D vortex patch problem, the velocity of a particle on the contour of 
the patch is obtained, following Zabusky et al. [Tl], by inverting the relation 
between the streamfunction and the vorticity. For the a-patches problem this 
relation is generalised hy 9 = (— A)-'^~^Vj where 9 is the scalar and is the 
streamfunction. In the limit a — > 0, which corresponds to the 2D Euler equation, 
the scalar plays the role of vorticity and in the limit a — ^ 1 (the quasi-geostrophic 
equation) the scalar corresponds to the potential temperature. For < a < 1 
one gets: 

«(^(7, t),t)^^ f -.-izlli!) (2) 



2^ Jc(t) K(7,0-^(7',i)l" 

where Izf (7, i) is the position of points over the contour C{t), parametrised with 
7, and 9i = 9ca. Here 9 is the value of scalar in the patch and the factor 

r(^) 

Cq, = — results from inverting the operator (— A)^"^. The contour 

dynamics equation is then obtained by stating: 

%^=.(^(7,t),0. 
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In this article we will focus on the contour evolution of two patches for which 
the equation become: 

dt ^ 27ric.(t) 1^(7,0 -^fe(7',t)l"' ^' 

Local existence for the limit a = 1 has been reported in [T31 [IS] . In this case 
one should use (see [13]) the following formula for the velocity (in the one patch 
case): 

«(^(7,t),i) = ^ / X J Jr J d^ (4) 



2^ Jc{t) \x{l,t)-lt{i,t)\- 

This equation eliminates the tangential component of the velocity, as only the 
normal component is able to deform the curve, and avoids therefore divergent 
integrals. 

Simulations of Eq. ^ as reported in show evidence of self-similar col- 
lapse. The singularity is point-like and we denote the collapse coordinates as 
!?*(<*), where is the collapse time. Numerical results of p] indicate that near 
the blowup there exist scaling laws in the evolution of the maximum curvature 
as 

C 



and also in the minimum distance between contours 

d^C{t^,~t)i as t ~* (5) 
We re-scale the spatial variable as 

^{j,t)-lt4t) = {t*-t)\^{j,t)^^,), (6) 
where 5 — 1/a, and the velocity of the scaling function ~a?,(i) is given by 

dic^t) Ok r 0il\t)dj' 

The equation of motion for the rescaled space variable if is: 

V^'Jt.c) Vl""!^-')- "■'(^'•'•I" llT.-fth'.*)!-/ 

Where are the rescaled contours. The new time variable 

r = -log(t, -t), (8) 
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transforms the equation ([7]) into the self-similar equation: 



k 2^ J^k(r) \\yiii'r)- y fc(7',T)l" \y*- yk{'y',T)\ 




(9) 



This equation, as the non-rescaled one, admits a projection over ~Un, the unitary 
normal component to the curve at 'y', which is useful to describe the a = 1 case, 



= iKy - 2/*))ir,. + 




, ^ 27r 



The projection evolves in time the same as the system ([9]) as only the normal 
component deforms the curve. Near the collapse time, rescaled contours re- 
ported in the work by Cordoba et al. (see [2j) seem to coincide over a unique 
curve, thus suggesting that Eq. (fTO|) has a fixed point, also called self-similar 
solution. 

There exist interesting relations between the original and rescaled variables. 
For instance, the blowup in the original variables is transformed into an asymp- 
totic behaviour in the self-similar variables. In this way the collapse reported in 
[2] that takes place in a very small time interval in original variables, is mapped 
to an infinite interval in rescaled variables, thus permitting a detailed moni- 
toring of the blowup. We consider now the area of the patches. It is related 
to the energy of the initial data and it is conserved in the original variables 
{Ic, t), however following Eq. ([6]), it tends to grow in the rescaled ("y*, r). Other 
connections between the rescaled and non rescaled variables regard to distance. 
The shortest distance d, between collapsing contours (see Eq. (O) tends to zero 
in the original variables but it is not necessarily zero in the new variables. In 
fact in [2] it is reported how the rescaled curves seem to keep a finite distance be- 
tween them. These results are easily justified. Since collapse is point-like, at the 
blowup time distance between contours becomes zero just for two trajectories, 
one on each contour, satisfying, 

\^cAt)-^cM-^o, t-^t* (11) 

Distance between these trajectories in self-similar variables, is given by, 

l^cAt) - ^cAt)\ = (t* - tYltcAt) ~ ^cM (12) 

The factor (i, — ty in Eq. confirms that the null distance between the 

original variables is satisfied even if trajectories on the rescaled variables are at 
a finite distance at the collapse time . The above expression may be rewritten 
in terms of the pseudo-time r as follows, 

I'^cAt) - ^cM = e-^'\tcAr) " yc.(r)| (13) 
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As explained in next section, the numerical technique used for simulations does 
not track individual particles on contours, but just the contour as a whole, so 
that the above collapse trajectories are not numerically integrated. However it 
is expected that near the collapse, the distance between these trajectories will 
be well represented by the evolution of the shortest distance between contours 
for which in [2] it is reported 

\itcAt)-^cAt)\-D-iu-ty/'', t~.U (14) 

Here lci,.(t),i = 1,2 represents the set of points on the contours that are at a 
shortest distance, as a function of time. Similarly in the self-similar variables 
the collapse may be tracked by means of the shortest distance between contours. 



e 



'yc^-lJc2\. (15) 



This expression indicates that existence of collapse in the self-similar variables, 
implies that the minimum distance between contours defined as: 

1^,, -^,J^^(r), (16) 

is a positive function, D{t), satisfying: 

hm ^ -> 0, (17) 

Otherwise the initial patches would be non collapsing. In other words, non 
collapsing patches, in the original variables remain always with finite curvature 
and at a finite distance at any time, but in the self-similar variables its minimum 
distance asymptotically goes at least as e'^^ , so that they are repelled one from 
other. Typically the time evolution of non collapsing initial data is not described 
in the self-similar variables framework, however this is done in this article guided 
by the exact self-similar solutions of the system for which numerics shows 
evidence that they constitute a boundary between collapsing and non collapsing 
data. We will come to these issues in detail in the results section. The existence 
of a T dependent function D{t), is consistent with the collapse classifications 
in terms of self-similar variables reported in |15| . There it is shown that the 
asymptotic behaviour of the collapsing data may be towards a fixed point (a 
constant D{t)) or towards r-dependent solutions which may be either periodic 
or chaotic. In this article we will show that numerical simulations indicate that 
the collapse of the a-patches problem address an asymptotic chaotic regime, 
suggesting that this could be the type of collapse present in the 3D Euler problem 
as conjectured in [TB]. 

3 The numerical method 

We have numerically solved Eq. ([3]) for < a < 1. The evolution is calcu- 
lated using ideas developed in 8, 9J for the case of the 2D Euler equations. In 
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X, (P) 




Figure 1: A representation of the contour Ck by a set of Nk points. The sohd 
Hne represents the interpolating curve between nodes x j and x j+i. 

particular, as described in each contour Ck is represented by a set of nodes 
Nk- Fig. 1 depicts a scheme showing how the curve is interpolated between 
consecutive nodes with a cubic spline, as follows: 

~x j {p) ^^^+p~tj+ (p) lij (18) 

for < p < 1 with ~a?j(0) — ~Xj and ~a?j(l) = where: 

7, = {aj,bj) = ^j+i-^j, eR2 (19) 

It J = {-bj,aj), e (20) 

Vjip) = ^ijP + PjP^ +"/jP^ Vj (21) 
The cubic interpolation coefficients fij, (3j and are: 

Mi = ^'^djKj — -djKj+i, (3j = 2djKj, — grfjlKj+i — Kj), 
where dj — |~a?j+i — ~Xj\ and 



(22) 



is the local curvature defined by a circle through the three points, 
and Xj+i. The node spacing in each contour is non-locally adjusted at each 
time step depending on curvature. Issues related to the density of nodes along 
the curve will be addressed later. Now it is explained how the system ([3]) is 
transformed by the above discretisation. The evolution of any point ~x{'-y,t) on 
the contours is replaced by the evolution of a point ~x j on the discrete curves, 
and the integrals on the contour curves Ck are replaced by the summation of 
integrals over the parameter p: 

dxjjt) ^ ^ / dp .231 
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here 'xi^k{p,t) refers to a piece of curve computed as in Eq. ([T5)) at a time t, 
and the additional subindex k distinguishes the contour where the segment is 
placed. More abstractly this expression may be written as: 



dt -'^^^^ 



f,{-^,{t),Ci,C2), /, :M2^R2 (24) 



This equation stands for any node ~x j on the discretised contour curves so, for 
the case of two contours on a plane, the system represents a set of M = 
2 X {Ni + N2) coupled ordinary differential equations, which more compactly is 
rewritten as an autonomous system as follows: 

^ = f(x), xGM^^ (25) 

The dimension M of the system (P5|) is very large and typically non constant, 
as the discretisation does not force to maintain the number of points Nk on the 
curve for all times. Following [9J the system (p5|) is integrated with an explicit 
^th Qj-^Qj. Runge-Kutta method. The time step At in the Runge-Kutta method 
as suggested by Cordoba et al. [5] is chosen dynamically, since near the collapse 
time it has to be refined. As it is reported there, the blowup concurs with the 
formation of corners on the contours, and the node spacing Ax is reduced to 
represent corners properly. The time step is then adjusted as At = BAx, where 
B is empirically tuned for different values of a. 

In the numerical simulation of the a-patches problem the most challenging 
part is the evaluation of each fj in Eq. (f24|) as required by the Runge-Kutta 
method at each time step. In order to explain how this is achieved, next focus 
is on the contour integrals along each cubic contour segment (see Eq. (|23)) V 
In this equation the interpolation is replaced and integrals are leaved as 
follows: 

\{t i+ ^iiUi) + (2/3ip + 37,;p^) ) dp 

^ ,^ ^ ^ , = (26) 

\x i- X J +pt .i + rii{p)ni\°' 

,^ ^ ^,f' dp _ mp + Sj,p-')dp 

[t 1+ fit Hi) / — — — ; — h n 



\x i - X J + p t i + 'i]i{p) Uil'^ Jo \x i - X j + p t t + rii{p) nil 

The subindexes k appearing in integrals of Eq. (j23|) have been dropped to 
avoid cumbersome notation without lost of generality. Attending to the relative 
position of the variable ~Xj versus the fixed ~Xi we classify the kind of integrals 
to do. 

1. The In this situation the first integral on the right hand 

side of Eq. ((26|) has an integrand that becomes infinity at p = 0. However in 
spite of this the principal value of the integral is defined if a < 1. In order to 
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compute the PV we rewrite the integral as follows: 



dp 



\p t i + r]i{p) n 
p~°'dp 



V + ) 

1 

p-" (C0 + C1P + C2/) 

. i V - (27) 

Coefficients on the expansion may be easily computed with Maple or Mathe- 
matica or other algebraic manipulator. In our numerical computations we have 
used an expansion up to the tenth order, which is very accurate as coefficients 
Cn rapidly decay. We provide an explicit expression for cq and ci . 



Co = 1, Cl = 



1 + 



These coefficients are also useful to approach the second integral on the right 
hand side of Eq. p6p . In particular this integral is approached by the summa- 
tion: 

|t,h(l + A^f)"/%eo "U-« + 2 + n-a + 3j' ^''^ 

2. The case x j — x ^+i. In this situation the first integral on the right 
hand side of Eq. (|26|) has an integrand that becomes infinity at p = I, as 
Xi{p = 1) = .Ti+i. Results described in 1. are applicable here if integrals are 
rewritten for p' = 1 - p, /I^ = + 2/3i + 87^, /3i = -fii - 871 and 7^ = 74. 

3. The case d^ = \~x j — ~x i\ » 0. In this case the first integral on the right 
hand side of Eq. (|26p may be approached by a series expansion as follows: 

dp 



\x r- X J +pt ^ + rii{p)ni\°' 
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Jo- / , ^ r, 19 ^ ^ '! 3" (29) 

where d„ = ^("^j ^ ~2?i) ■ and dt — — ~x'i) ■ t , and coefficients are 

computed with a software for algebraic manipulation. Coefficients g„ rapidly 
decay if the distance, dx, between Icj and Ici is large enough. We provide an 
explicit expression for go and gi: 

_ , _ a{dt + dnUi) 
2d? 
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Figure 2: A graphic of the integrand in the first integral at the right hand side 
of Eq. (UHl) in the 1^0. It may vary steeply in the domain. 



The second integral on the right hand side of Eq. ([26| may also be expressed 
in terms of coefficients (/„. 

4. The case I'j^j — ~ 0. In this situation the first integral at the right 
hand side of Eq. has an integrand w{p) which for instance may take the 
appearance depicted in Fig. 2. As noted there, it may vary steeply in the 
domain. In this case we compute the area below the function: 



/ w{p)dp — A 
Jq 



by integrating the following ordinary differential equation: 

dY 

— = w{p) 
dp 

with initial condition Y{p = 0) = 0. With this choice the integrated function 
will be such that Y{p = 1) = ^4. The integration method used is a variable step 
e^th Qj-jgj. Runge-Kutta method (see Tf]), thus ensuring it will appropriately 
track steep variations. The evaluation of wij)) must be done with a lot of care 
to avoid round off errors. This is correctly done if the denominator in Eq. (|26p 
is computed as follows. First the polynomial in p is expanded and the coefficient 
of each monomial in p is evaluated. Second these terms are added up. 

In the numerical code an appropriate test distinguishes which one of these 
procedures needs to be followed to evaluate an integral. The code is verified by 
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evaluating fj{^j, Ci) for a case with only one contour, the circle of radius 1 cen- 
tred at the origin, afj = (1,0) and 9i = — 1. The output of the discrete scheme 
just explained is compared with that obtained from the natural parametrisation 
in 0, i.e. Ci = (cos(6'), sin(0)), 9 G [Q,2tt]. In this case the integral to evaluate 
is: 

-1 f^"" (-sin 6*, cos 6*) de 



^/(l - cos 0)2 



sin^ ( 



This integral is numerically evaluated with Maple forcing 5 digits precision, 
which gives fj('xj,Ci) = (9.549296 • 10-l^ -8.4000655 • lO^^). The above 
numerical method considering 200 points on the circle supplies: fj{licj,Ci) = 
(-3.09849665 • 10-^ -8.39990851 • 10"^). Both results agree within reasonable 
accuracy. The precision of the numerical method depends on the number of 
points on the discretised curve, and also on certain tolerance in the criteria used 
to distinguish between the application of the procedure 3 or 4 discussed above. 
In practise evaluation of fj is accurate up to the 5*'' digit in consistency with 
the 4*^* order Runge-Kutta method employed to evolve x. 

The evolution of the quasi-geostrophic equation (i.e the case a = 1) as 
reported in Eq. ([4]) is computed projecting the vector field of Eq. ([3]) over the 
unitary normal component to the curve at ~Xj. This means that numerically 
the first integral at the right hand side of Eq. ([26]) needs not to be computed 
(in fact it is a divergent integral (see [T^)) as it is a tangential component to 
the curve at x j. The second integral at the right hand side of Eq. (|26p need to 
be projected over the unitary normal component to the curve, and this is easily 
done. The evolution of the contours for a < 1 may be indistinctly computed 
either with Eq. ^ or with Eq. (g]). 

Once f (x) in Eq. (|25p is conveniently approached we are ready to apply the 
/^th Qj.(jgj. Runge-Kutta method to step forward x. Following the methodology 
developed in [9] every time step the nodes are redistributed on the contour to 
guarantee its optimal representation. This means that the numerical method 
will not track individual trajectories, but the whole contour. The relocation of 
points is done as follows. For each node assume a local density 



1 + <5k,/V2 

This form limits the minimum distance between nodes to nearly 5 and the 
maximum computed curvature to . We take 5 = 10~^. Here ki = [hi -f 
Kj+i)/2 and 

ki = {vLy^{k^LY + ^/2k^ 

We have considered L = 3, a = 2/3, and v \s a. parameter which controls the 
accuracy of the representation, being 0.05 or 0.03 typical values for v. Smaller 
values imply larger densities of points on the curve. The curvature k incorpo- 
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rates non-local effects and is defined as: 




where dj = — 'x'j\, hij = i — ("a^j+i + l?j)/2| and Rj — {kj + Kj+i)/2 

(see Kj defined in Eq. (f22|) ) and the summation is over all nodes in the contour 
of the node i. 

Given a desired density pi between nodes i and i + it is possible to compute 
(Ti = pidi, which is the fractional number of nodes to be placed between i and 
i + 1. The node redistribution in a contour k fixes one point on the curve, 
the first point, which will be kept fixed on the redistribution procedure. This 
implies that this point will evolve as a real trajectory of the system ([25]) . The 
quantity q is computed: 

and Nk=[q\ + 2, the nearest integer to q plus two. The Nk — 1 'old' nodes 
on the curve, will be replaced by Nk — 1 new nodes which will lie along the 
curved contour segments connecting the old nodes. Let cr- = aiNk/q, so that 
^i^i'^'i = ^k- Then the positions of the new nodes j — 2...Nk are found 
successively by seeking i and p such that: 

i-i 

Y.<y'i + <y',P = j^l (30) 

and placing each new node j between the old nodes i and z + 1 at the position 
"x i(p) given by Eq. (fTS]) . 

The evolution of the self-similar equation ([9]) may also be computed numer- 
ically, as its discrete version may formally be rewritten as: 

i _ 17. . TP. ■ 1D)2 



dr 



= F,(^„Ti,T2), :R^^R^ (31) 



The numerical evaluation of Fj at each step of the 4*'* order Runge-Kutta 
method may be easily achieved with the use of the algorithm developed for 
the fimction fj as it is satisfied that: 

F,(y„Ti,T2) = <5y, +/,(^„Ti,T2) - /,(0',Ti,T2) 

where the choice if^ = has been considered. 

4 The results 

The numerical scheme explained above has been applied in [2] to simulate the 
time evolution of two identical circular patches of unit radius with centres placed 
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at a distance 2.5 and such that 0^ = — 1 for /s = 1, 2. In particular Eq. ^ is 
studied for 0.5 < a < 1 and Eq. ([4]) for the case a — 1. In Ref. [2] it is provided 
evidence of the evolution of the initial data towards a collapse presenting a rather 
stationary rescaled profile. Simulations performed in [5] were not implemented 
in the self-similar variables, and as novelty the results discussed in this article 
adopt the perspective of the self-similar variables. This section discusses the 
exact and numerical solutions to Eq. ([9]). Exact results are valid in the range 
< a < 1, which encompasses the quasi-geostrophic equation. Numerical 
simulations have been just obtained for a ~ 0.7, although they are expected to 
reflect the evolution on a whole a-range. 



4.1 The self-similar solution 

There exist exact self-similar solutions for contours Ti and T2 in Eq. (fTO|) . As 
expected from rescaling (jG]) at the collapse time the area within self-similar 
contours is infinite. Stationary-like profiles of the rescaled curves presented in 
[2 suggest that each of these curves take the form of a function T — {x,y{x)). 
When curves are parametrised as functions Eq. (fTU]) is rewritten as: 



Ot ^ ^ 27r 

^ix',T)dx' I ^{x',T)dx' 



I 



^ix,r)dx' ^ r ^{x',T)dx' , 
'y*-^i{x',T)\" Jr IV* - V2(a;',r)|«- 



(32) 



A fixed point to Eq. is: 

Ti = (x,|x|), T2^{x,-\x\) (33) 

These curves touch each other at one point, what is compatible with results on 
regularity of the contours discussed in [18] . Figure 3 shows geometrically how 
over these particular curves, the contributions of different terms of Eq. (|32p 
to a point, marked with a black dot, cancel each other out. The first term of 
Eq. p2p evaluated at the black dot, provides a null contribution on the normal 
direction, as the arrow labelled with '1', is tangent to the curves ([33|) . Arrows 
labelled with '2' and '3' show the normal component of the small contributions 
of the second and third terms of Eq. ([5^ integrated along the small segments 
within the gray dots. They cancel each other as they are equal magnitude but 
opposite sign. Similarly contributions from the 4*'* and 5*^ terms cancel each 
other. 

The fixed point in Eq. (|33p is not unique, as any rotation of it around the 
origin is also a stationary solution of Eq. (jlOp . There exist infinite stationary 
solutions, related each other by this continuous transformation. There exist 
previous examples in the context infinite dimensional dynamical systems, of 
families of fixed points related by means of a group action (see for instance 
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Figure 3: Sketch showing contributions of different terms in Eq. (|32p to the 
point distinguished with a black dot. Each arrow is labeUed with the number 
of the term it corresponds. It is clear that the geometry of the curves imposes 
cancellations of the normal components. 

[20]). It has been reported [20J that this transformation implies the existence of 
a null eigenvalue along a dynamically neutral eigenspace. Apart from these fixed 
points, there might be other stationary solutions to Eq. IjlOp . however we have 
not found those. Fixed points in the theory of dynamical system are particularly 
useful solutions, as they geometrically organise the flow on the phase space. Eq. 
(jlOp can be viewed as a dynamical system of infinite dimension that can be 
represented by a large number of equations as in Eq. (pS)) . The fixed point 
in Eq. ([55)1 is not a global attractor of pI7|) . This could not be the case as 
otherwise every initial condition would be attracted by it, and this would imply 
non existence of non collapsing solutions. On the hand conservative systems, as 
is the case of the non rescaled system, do not have solutions which are attractors. 
Numerical simulations confirm the existence of perturbations to the curves (|33p 
which are completely unstable. Figure 4a) shows the backwards time evolution 
of an initial condition which is a small perturbation to the self-similar functions 
(|33|) . This initial data after 14 time steps with negative At approaches more and 
more the fixed point, implying that this perturbation is completely contained in 
the unstable subspace of . Figure 4b) displays the evaluation of the modulus 
of the normal component of the velocity along these curves and as expected it 
is almost zero. In these simulation curves Ti, T2 have been transformed in 
functions with the asymptotic behaviour at large x extended with straight lines 
following the self-similar behaviour. These lines are large but not infinity, thus 
errors are introduced at their boundaries, however these errors are very small 
in the displayed zones. 
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Figure 4: a) Initial data (solid line) and its backwards time evolution (line with 
circles) according to Eq. pop : b) Modulus of the normal component of the 
velocity along the initial data curves. Solid thick line corresponds to the lower 
curve, while thin line corresponds to the upper curve. 
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y, y. 



Figure 5: a) Initial data and the self-similar solution; b) Repelling r evolution 
of the patches in a). 

4.2 The separatrix and the 'chaotic' collapse 

The evolution of curves Ti, T2 is described in self-similar variables. As it has 
been already pointed out, from Eq. ^ follows that in these coordinates the 
area enclosed by the contours becomes infinite at the collapse time t,, but on 
the other hand blowup takes place at infinite time in the pseudo-time t, as 
transformation ([8]) confirms. This implies that at any stage of the simulation of 
Eq. ([9]), T and also the area of the patches are both finite, making possible the 
numerical computation of its evolution. 

The solution (|33p divides the plane in four domains that are labelled as the 
upper (U), the bottom (B), the right (R) and the left (L). The distribution of the 
patches over these regions provide configurations of collapse. In fact numerical 
evidence is found on the role of the fixed point ([33]) acting as a separatrix 
between non collapsing and collapsing initial data. In particular as discussed 
next a non collapse configuration of the initial data cannot cross the boundaries 
of the fixed point. 

Figure 5a) displays two patches that initially stay very close each other, but 
at the same time they closely follow the fixed point ([33]) plus a small perturbation 
in the unstable subspace. This initial data is entirely contained in the upper 
(U) and lower (L) domains. Figure 5b) confirms the repelling r evolution of this 
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Figure 6: r evolution of the logarithm of the minimum distance between rescaled 
contours. Crosses stand for non coUapsing initial data. The dashed linear 
fitting of the last N data provides a gradient m = 1.4742 which confirms a non 
collapsing condition according to the limit ^T7} . The black and gray solid lines 
respectively correspond to the evolution of the collapsing initial data of Fig 8a) 
and 8b). 

initial data. The minimum distance, D{t), between these patches evolves with 
T as indicated by crosses in Fig. 6. The logarithmic scale on the vertical axis 
confirms that the r evolution is in close correspondence to Eq. (fT3)) . A linear 
regression adjusting = 16 data in the interval t e [1.5, 5.5] provides a gradient 
m = 1.4742 very close to the limit of the collapse condition S = 1/a — 1.4286 
provided in (jl7p . In particular non collapse is guaranteed if m > S. Similar 
results would have been expected for any arbitrary rotation of the initial patches. 

Simulations reported in j2] address the existence of collapsing initial data. 
In particular there it is reported that two circles in which 6k = —1, with unit 
radius and with centres placed at a distance 2.5, collapse in a finite time. The 
collapse is point-like and the symmetry of the initial data reports two collapsing 
points on each curve. Blow up also exists for initial data different from that 
discussed in [5] . For instance if the above described circles are deformed to two 
ellipses with semi- axis a = 1.1 on the horizontal line and 6 = 1 in the vertical 
line, the collapse regime described in [5] is also reached. The collapse time of 
these ellipses is estimated by means of a nonlinear adjustment as explained in 
[2l[T9] and it is obtained that in this case = 6.887794662. This value is used to 
rescale profiles according to expression ([6]) . In particular we perform simulations 
in the self-similar variables taking as initial data the curves of this simulation 
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Figure 7: r evolution of the logarithm of the minimum distance between rescaled 
collapsing contours. The black line correspond to a simulation with low reso- 
lution V — 0.05 and the gray line corresponds to the same initial data with a 
higher resolution v =^ 0.03. Both resolutions confirm a highly oscillating system. 

at time tmax = 6.81482306, and rescaling them with the transformation In 
this expression y = (0, 0) and is an unknown rescaling function. Figure 

8 shows the two choices taken for it. In figure 8a) it is considered ~x *{tmax) = 
(0,0) and in figure 8b) it is considered Ic ^{t,nax) = (1.36439899,-0.28345455) 
that moves the two collapse corners to the left (L). In Fig. 8a) the configuration 
is as follows: one patch overlaps regions (L), (U) and (R) while the other overlaps 
regions (L), (D) and (R). For the Fig. 8b) the configuration is: one patch is 
contained in region (L), while the other overlaps regions (L), (D) and (R). In 
both cases the minimum distance between patches evolves according to the 
collapse condition (fT7| . The function log(Z)(r)) is depicted in Fig 6. The 
solid black line stands for the initial data of figure 8a) while the solid gray line 
corresponds to the initial data in figure 8b). Both regimes stay below the line 
which mark the limit for non collapsing initial data. Simulation are performed 
for a T interval in which the area of the patches increases 4 orders of magnitude. 

Figure 7 shows a magnification of the evolution of log(-D(r)) for the initial 
data in Fig. 8a). It is clear the rapidly oscillating behaviour that does not 
settle in a stationary regime as would be the case of an attracting fixed point. 
The black line corresponds to the evolution of the minimum distance with res- 
olution V = 0.05, the gray line corresponds to a higher resolution (y = 0.03). 
In the latter case the amplitude of the oscillations diminishes, but still con- 
firms a "chaotic" evolution of the distance. It seems that these oscillations are 
under-resolved in both discretisations, and there is a noticeable sensitivity in 
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Figure 8: a) Initial collapsing data; b) A second initial collapsing data. 
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Figure 9: Initial data which collapse slowly. 

the evolution to the dimension of the discrete system ([25l) which is different in 
both cases. This confirms the conjecture of a "chaotic" evolution embedded in 
an infinity dimension dynamical system. 

Fig. 9 shows another collapsing data distributed like in Fig. 8a). However 
in this case the non collapse configuration discussed in Fig. 5 {i.e. each patch is 
entirely contained in opposite domains) is only slightly violated. Numerically is 
found that the evolution of the distance between contours follows an exponential 
law as the initial data of Fig. 5 did. However now the gradient obtained from 
the linear fitting t vs log(I?(T)) is to = 1.4187. This is almost equal to S, 
however it does not satisfy the non collapse requirement m > 6, and thus the 
initial condition provides an example of a slowly collapsing data. 
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5 Conclusions 



We have studied the a-patches problem in self-similar variables and we have 
found a family of exact self-similar fixed points for these equations. These are 
valid for the whole parameter range < a < 1 which encompasses the surface 
quasi-geostrophic equation. For the a = 0.7 case intensive numerical simula- 
tions are performed. They indicate that the stationary solution organises the 
evolution of different initial data, since it acts as a separatrix of collapsing and 
non collapsing conditions. It is found a perturbation to this solution that com- 
pletely fall in its unstable subspace. These are repelling patches that separate 
exponentially at a rate m > S thus confirming the non collapse requirement 
(fT7|) . The T evolution of other initial conditions is also discussed. There are 
found collapsing solutions for initial configurations that cross the boundaries 
marked by the fixed point ((33| . It seems that some of these data do not reach 
a stationary regime but a "chaotic" one. Also initial conditions that slowly 
collapse are reported. 
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